function [gbest,gbestval,fitcount] = QGSA_func(fhd,Dimension,Particle_Number,Max_Gen,lowerbound,upperbound,pos,bool_Plot,varargin)

g = 0.013;
rand('state',sum(100*clock));
me = Max_Gen;
ps = Particle_Number;
D = Dimension;

[pub, plb, vub, vlb] = SetBound(D,ps,lowerbound,upperbound);
% pos = VRmin+(VRmax-VRmin).*rand(ps,D);

fitcount = ps;

fitness = feval(fhd,pos',varargin{:});
sort_fitness = [fitness' (1:ps)'];
sort_fitness = sortrows(sort_fitness, 1);
fitbest = sort_fitness(1,1);
fitworst = sort_fitness(ps,1);
gbestindex = sort_fitness(1,2);
gbest = pos(gbestindex,:);
gbestval = fitbest;    

if(bool_Plot)
    figure(4)
    haxes = plot( 0 , 0 );
    XArray = [1];
    YArray = [gbestval];
    title('QGSA');
end

for iter = 2:me
    if(fitbest == fitworst)
        mass = repmat(1/ps,ps,1);
    else
        mass = (fitness'-fitworst)./(fitbest-fitworst);
%         mass = zeros(ps,1);
%         for n = 1 : ps
%             mass(n,1) = (fitness(n)-fitworst)/(fitbest-fitworst);
%         end
    end
    mass = mass / sum(mass);
    sort_mass = [mass (1:ps)'];
    sort_mass = sortrows(sort_mass, 1);

    numKbest = ceil(ps * (1 - (iter /me)));
    numKbest = (numKbest>1)*numKbest + (numKbest<=1)*1;
    pos_Kbest = zeros(numKbest, D);
    probability_Kbest = zeros(numKbest, 1);
    for k = 1 : numKbest
        index_Kbest = sort_mass(ps-k+1, 2);
        pos_Kbest(k,:) = pos(index_Kbest,:);
        probability_Kbest(k,1) = mass(index_Kbest);
    end
    
%     index_center_map = RouletteWheel(probability_Kbest, ps);
    for n = 1 : ps
        index_center = RouletteWheel(probability_Kbest, 1);
        center = pos_Kbest(index_center,:);
        for d = 1 : D
            S = rand();
            S = (S>=0.5)-(S<0.5); % S = (S>=0.5)*1+(S<0.5)*(-1)
            pos(n,d) = center(1,d) + S * g * abs(center(1,d) - pos(n,d)) * log(1/rand());
        end
    end
    pos = ((pos>=plb)&(pos<=pub)).*pos+(pos<plb).*plb+(pos>pub).*pub;
    
    fitness = feval(fhd,pos',varargin{:});
    sort_fitness = [fitness' (1:ps)'];
    sort_fitness = sortrows(sort_fitness, 1);
    fitbest = sort_fitness(1,1);
    fitworst = sort_fitness(ps,1);
    gbestindex = sort_fitness(1,2);
    gbest = pos(gbestindex,:);
    gbestval = fitbest;
    
    if(bool_Plot)
        XArray = [ XArray iter]; 
        YArray = [ YArray gbestval];
        set( haxes , 'XData' , XArray , 'YData' , YArray );
        drawnow
    end
end

%     if(mod(iter, 10) == 0)
%         iter
%     end
end
